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1. Introduction 

Frequency conversion, via three wave mixing (TWM) pro- 
cesses in quadratic nonlinear optical media, is widely used 
in order to generate laser frequencies that are not available 
by direct laser action ^ The efficiency of a TWM process 
depends on the fulfillment of a phase-matching condition''^. 
Quasi-phase-matching (QPM)'"^, a method in which the sign 
of the nonlinear coefficient is modulated, facilitates control 
over phase-matching conditions. Still, QPM processes are 
generally not robust against variation in system parameters, 
such as temperature, input wavelength, incidence angle, etc. 

Recently, several works have been published that concern 
robust adiabatic TWM processes in the fixed (undepleted) 
pump approximation''"'^, i.e. when one of the waves is much 
more intense than the others, and thus is negligibly affected 
by the interaction. This assumption linearizes the dynamics, 
making it isomorphous to the linear Schrodinger equation of 
quantum mechanics, and thus allows the use of quantum me- 
chanical adiabatic theorem'^. 

The first step towards fully nonlinear TWM was taken by 
Baranova et al.'"*, for the special case of second harmonic 
generation (SHG). Phillips et al. extended the work into the 
realm of optical parametric amplification (OPA) and optical 
parametric oscillation (OPO)'^!'^^. However, these works do 
not provide a rigorous physical model explaining the observed 
phenomena. Rather, it was stated that this is a generalization 
of the case with fixed pump, analogous with a quantum model 
of a two-level atom'^. This generalization is not self-evident, 
as the removal of the fixed pump approximation invalidates 
the analogy made with other systems. Specifically, a reference 
was made to the geometrical representation of TWM made by 
Luther et al.'^ as being analogous to that made by Crisp- 
with regards to a nonlinear two-level atom, which builds on 
the Feynman, Vernon and Hellwarth model We maintain 
that this analogy does not hold, since the nonlinearities in the 
two physical systems, TWM and two-level atom, are of differ- 
ent nature. The dynamics of the two-level atom remains linear 
at all times, as the effective wave vector is governed entirely 
by the electric field, which is taken to be independent of the 
atomic state in the approximation made by Crisp. The nonlin- 



earity is expressed in the resulting susceptibility of the atom. 
Contrarily, in the TWM geometrical representation, the anal- 
ogous quantity to the effective wave vector is a function of the 
interacting field amplitudes, which renders the dynamics itself 
nonlinear. Two exceptions are special cases for which a sound 
physical model was found: (i) the case studied by Longhi^*^, in 
which SHG was followed by sum frequency generation (SFG) 
to generate the third harmonic, which was found to be analo- 
gous to a certain nonlinear quantum system^' (ii) the case of 
OPA with high initial pump-to-signal ratio, which Yaakobi et 
al.^^ approached as a case of auto-resonance. 

Other groups have taken up quantum systems with fully 
nonlinear dynamics, and developed a theory of adiabatic evo- 
lution for them^^"^^ . Interestingly, they base their method on 
representing the Schrodinger equation in a canonical Hamil- 
tonian structure, as was done in classical mechanics, and use 
classical adiabatic invariance theorem"^. The equations gov- 
erning TWM have also been put in a canonical Hamiltonian 
structure in several worksJ^i^iS^, but not in the context of adi- 
abatic evolution. 

Here, a comprehensive physical model of fully nonlinear 
adiabatic TWM is presented for the first time to the best of 
our knowledge. This analysis leads to a condition for efficient, 
broadband and robust frequency conversion. Such conversion 
is demonstrated numerically. 

This paper is organized as follows. In Section |2] the theo- 
retical model of TWM is presented, along with this system's 
stationary states, using canonical Hamiltonian structures. In 
section III adiabatic evolution is analyzed, and an analysis of 
robustness leading to large bandwidth is provided. Section |4] 
presents numerical simulations of adiabatic TWM with phys- 
ically realistic parameters, available with current technology. 



2. Theoretical Model 

A. Coupled Wave Equations in Canonical Hamiltonian 
Structure 

The dynamics of TWM is commonly described by three cou- 
pled wave equations. Assuming plane-waves and a slowly 
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where yj — {kjC^^ are the coupHng coefficients, and 

kj and Aj are the wavenumber and complex amplitude of the 
wave at frequency 0)^, respectively. is the second order 
nonlinear susceptibility and Ak = k\ + k2 ~ ko, is the phase- 
mismatch. Without loss of generality we assume Oi < 02 < 
o>i where 0)3 = ©i + a>2. 

From this point, we follow the analysis of Luther et 
al3 in the construction of a canonical Hamiltonian form 
of the coupled wave equations. First, we define Aj ~ 

y/Yjqjexp (^—i jQAk(z')dz') and note that this renders \qj\^ 
proportional to the photon flux at coj. Next, we write the three 
equations using qj, 



d^ 
dqi 
'd^ 

dqT, 
'd^ 



iAYqi — iq*2q3 
iATqi - iq\qi 
/Ar^3 — iq\q2 



(2) 



where we also defined the scaled propagation length ^ — 
Z^/Yxjiji and the parameter AF = Ak/ ^yxy^Ji, which de- 
scribes the relative strength of the phase-mismatch compared 
to the nonlinearity. The coupled equations can now be written 
in a canonical Hamiltonian structure. 



dqi 
'd^ 



.dH 



(3) 



where qj play the role of the generalized coordinates, q*j are 
their conjugate generalized momenta and 



H ^ ij^{q\q\qi+ qmql) 



AF^ 



(4) 
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is the Hamiltonian. Additionally, we have the Poisson brack- 
ets relations 



which can readily be shown to be generated by the Manley- 
Rowe relations. 



^^1 = \qi 
K2 
K3 



K2 = \qi\^-\q2\^ 



111 



to I 



(9) 



i.e. the Kj are constants of the motion. 



B. Stationary States 

The stationary states are very significant for the adiabatic evo- 
lution analyzed in section III. It will be shown there that when 
an adiabaticity condition is satisfied, the system evolves along 
these states as they follow a slowly changing system param- 
eter - the phase-mismatch. Determining the dependence of 
these states on phase-mismatch is thus crucial for predicting 
the outcome of adiabatic evolution. 

The TWM system is known to have two stationary states— 
besides the trivial ones, i.e. the states where two of the three 
waves have no energy. For completeness, they will be derived 
here as well. We note that any parametric instabilities are ig- 
nored here, as we seek only stable solutions. 

In a stationary state, the state of the system is transformed 
into itself by the evolution dynamics. The coupled wave equa- 
tions |2] are invariant with respect to the transformation 

[q\,q2,q3.) {q\exp{i9i^) ,q2exp{i92^) ,q3exp[i{9i + 62)^]) 

(10) 

as evident from Eq. |6]and[8] If 



dqj 
~d^ 



^ = + 02)^3 



(11) 



{qi,qj} - 0, {q*,1*j} - 0> {qi,1*j} - -^idij (5) gjjj^pjg solutions: 



then Eq. |2]will perform the transformation [TO] and remain in- 
variant, i.e. the system state will be transformed into itself. 
Therefore, Eq. [TT]define the stationary states for this system. 
Substituting these relations in Eq. |2] yields quartic equations 
of 01 and 02> with the Manley-Rowe relations as parameters. 
For any given pair of Manley-Rowe constants, there exist 9i 
and 02 that yield two nontrivial stationary states, which we 
hence term the "plus state" and "minus state", and use cor- 
responding indexes in mathematical expressions. These so- 
lutions are very involved algebraically, and do not facilitate 
physical insight. We therefore focus first on the special case 
where the two low frequencies have the same photon flux, i.e. 
'qi\^ = Iqil^ (note that still, generally, coi 7^ O2), which leads 



Finally, we note that the Hamiltonian is invariant under the 
phase transformations 



{q\-,q2,qi) {q\exp{i9\),q2,qiexp{i9i)) (6) 
{qi,q2,q3) ^ {q\exp{i92) ,q2exp{~i92) ,qj.) (7) 
{q\,q2,q?) {q\,q2exp{i9i) ,q^exp{i9i)) (8) 
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Fig. 1. Normalized photon flux of each wave of the two sta- 
tionary states with \q\^ = |^2p- (a) Minus state, (b) Plus 
state. 



Fig. 2. Normalized photon flux of each wave of the two sta- 
tionary states with K2/ {Ki +K-i) = 0.3. (a) Minus state, (b) 
Plus state. 



and 



?1 ^92 



f v/(Ar-0_)(Ar-20_)eji:p ^ ) , AF < y/IP^ 

[0 ,Ar ^ 

f (AF-0_)ex/?(2/0_^) ,Ar<V2i^ 



where 



0± 

^3 



(14) 



The normalized photon flux of each of the three waves, as a 
function of the normalized (dimensionless) phase-mismatch 
AT/y/Pi, for each of the stationary states, is plotted in Fig. 
[T] for the case where \qi\ — \q2\ ■ For the minus state, as 
Ar/^/p5 approaches —00, the photon flux of the waves with 
the two lower frequencies (i.e. (Oi and CO2) approaches P3/2. 
It monotonically decreases with increasing AF/^/P? up to 
AF/ y/P^ = V2, where it vanishes and stays nulled for any 
AT/y/P^ > y/2. The high frequency wave (03) photon flux ap- 
proaches for AF/ —00, monotonically increases with 
AF/V^ up to AF/\/^ = \/2, and stays constant at P^ /2 for 
any AF/P3 > y/2. The dependence of the plus state intensities 
on AF is the mirror image, around AF — 0, of the minus state's 



intensities dependence, i.e. 



that where \qf\^ 



m I 



Note 



the stationary states are in fact 
trivial. 

Fig. |2] shows the photon flux of each wave of the sta- 
tionary states with the same parameters, for the case where 
7^ l'?2p- For the minus state, the three waves have the 
same monotonic dependence on AF/ ^/Pi as in the special case 
of \qi I — \q2\ , except that the two low frequency waves do 
not vanish (the kink that was observed at AF/^/P^ = y/2 in 
Fig. [H is now missing). Instead, of these two waves, the 
one that has the lower photon flux ( | in Fig. |2^) asymp- 
totically approaches zero, while the other one remains at a 
constant difference from it, which corresponds to the value 



of K2 that characterizes this state. Since the stationary state 
is also characterized by a certain value of ^^3, ^3 always 

complements Iqjl to maintain the same ^^3. These sta- 
tionary states are thus never trivial. Furthermore, as before, 

^t(AF)p==|^7(-AF)|'. 

C. Dimensionally Reduced Canonical Hamiltonian Struc- 
ture 

The two previous subsections summarized representations 
and properties of the TWM system that were already known. 
In this subsection, a new representation is developed. This 
representation will be used in section III to account for adia- 
batic evolution. 

The existence of the constants of the motion Kj, in addi- 
tion to H, indicates that the number of degrees of freedom of 



the system is lower than the dimensionality of the [qjAj 

phase-space. As Liu et al.^^ have done for systems with t/ (1) 
symmetry, we'll use these constants to produce a phase space 
with reduced dimensionality. We define the real generalized 
coordinates Qj and real generalized momenta Py. 

Q\ = ~^arg{qi)-^arg{q2) + ^arg{qi) 

Ql = -^arg{qi) + ^arg{q2) 

Qi = -iarg{qi)-^arg{q2)-'^arg{q3) (15) 



Pi 

Pi 
P3 



|2 



\qi 

K2 = 1^1 

Ki+K^ = \qi 



|2 
121 

2 



-2 1^3 



(16) 



Ql is proportional to the phase difference between the two low 
frequencies and the high frequency, Q2 is proportional to the 
phase difference between the two low frequencies, and is 
proportional to the sum of phases of all three waves. Corre- 
spondingly, Pi represents the excess of photon flux in the two 
low frequency waves over the high frequency wave, P2 repre- 
sents the excess of photon flux at cOi over 02, and P3 repre- 
sents the overall photon flux balance between the three waves. 
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Using these definitions, the canonical Hamiltonian wave equa- 
tions become 



dQj _ dH dPj _ dH 



(17) 



with the Poisson relations 

{Qi,Qj} = 0, {Pi.Pj} ^ 0, {Qi.Pj} - d,j (18) 
and the Hamiltonian 

H = ^ViPi+ 2P2 + Pi) (Pi - 2P2 + Ps) i-Pi + Pi)cos (Sgi ) 

(19) 



-— (P1 + 3P3) 



The Hamiltonian is independent of Q2 and Qt,, indicating that 
P2 and P3 are constants of the motion, which is not surpris- 
ing since P2 = K2 and Pt, = K\ + Kt,. P\ and Qi thus form 
a closed set of Hamiltonian dynamics. We further note that 
the simple requirement that \qj\ > 0, 7 = 1,2,3 results in 
limiting the range of physically significant values of P\ to 
2 \P2\ -Pi<Pi< P3, for given P2 and P3. In fact, this ex- 
actly corresponds to the range of Pi for which H is real. Note 
also that, since Pi is bounded from below by 2 |P2 1 — when 
P2 7^ it sets a limit on the minimum value of Pi . This can 
be understood from a physical point of view: if P2 7^ then 
the photon fluxes at coi and a>2 are not the same. In upcon- 
version, each photon contributed to 0)3 by one of these waves 
is accompanied by a photon from the other wave, and causes 
Pi to decrease. When one of these waves is depleted the up- 
conversion process cannot continue, so Pi can no longer de- 
crease. When either \qi |^ = or \q2\^ = then, by definition. 
Pi = — 2P2 — P3 or Pi — 2P2 — P3, correspondingly. We further 
note that for P2 and any finite AT, Pi = 2 IP2I - P3 is a 
trivial stationary state, however for |Ar| 00 it is not. 

The stationary states correspond to fixed points in the 
(2i ,Pi ) phase space where 



dQi 

dPi 
1^ 



dH 



= 



dH 



The second equation results in 



= 



(20) 



(21) 



For the special case where the two low frequencies have the 
same photon flux 

P^ = 2|(Ar-0±)(Ar-20±)|-2(Ar-0±)2 (22) 

and the constants of motion P2 and P3 take the values 
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Fig. 3. Phase space portrait with normalized phase-mismatch 
AF = 0.6y/P^ and (a) P2 = (b) P2 = O.3P3. Arrows indicate 
motion of fixed points with increasing AF. 




Normalized Phase-Mismatch Normalized Phase-Mismatch 

Fig. 4. Normalized reduced phase-space canonical momen- 
tum for the two stationary states with (a) P2 = (b) P2 = O.3P3. 



Fig. [3^ and[3}5 show the reduced phase space portrait with 
P2 = and P2 O.3P3, respectively, where in both cases the 
phase-mismatch is AF = 0.6^/Pi. The fixed points, which 
correspond to the stationary states, are labeled by their in- 
dexes. The arrows indicate the direction of motion of the fixed 
points with increasing phase-mismatch AF. Fig. |4^ and ^ 
display Pf' as a function of the normalized phase-mismatch 



2|(AF-0±)(AF-20±)|+2(AF-0±)^ (23) 



AF/v^ for each of the two stationary states, with P2 = 
and P2 = O.3P3, respectively. Fig. shows that Pf P3 for 
AF ^ — and it decreases monotonically with increasing 
AF up to AF = \/2^. For any AF > \/2Pi, it stays constant 
at — P3, all in correspondence with the intensity dependence 
shown in Fig. [1] Similarly, Pj+ is the mirror image of P^ 
around AF = 0, i.e. P+ (AF) = Pf (-AF). In Fig. gb it is seen 
that Pj^ have the same monotonic dependence on AF/y^ as 
in the P2 = case, except that it persists throughout the en- 
tire range of AF/y^, i.e. there is no kink as in the previous 
case. Instead, with increasing AF/ \/p3, Pf goes from P3 to an 
asymptote approaching 2 IP2I — P3, and Pj is its mirror image, 
as before. 

3. Adiabatic Evolution and Bandwidth 

A. Adiabatic Evolution 

According to classical mechanical theory^^, an elliptic fixed 
point will follow an adiabatically varying control parameter, 
i.e. a parameter that changes slowly compared with the fre- 
quencies of periodic orbits around the fixed point. It will be 
shown how this adiabaticity condition naturally arises from a 



5 



linearization of the canonical Hamiltonian dynamics, i.e. Eq. 
[TTl about the fixed point^^'^^, where the adiabatically varying 
parameter is the phase-mismatch AT. The main result of this 
work is the derivation of the adiabaticity condition, as will be 
outlined below. 

The linearization procedure of Eq. [17] is detailed in ap- 
pendix A. It is shown that the nontrivial stationary states cor- 
respond to elliptic fixed points, and that 



5P, 



(24) 



where 5Pi = f i — f [ , i.e. it is the vertical difference between 
the system point and a fixed point in the {Qi,Pi) phase-space. 
V is the frequency of periodic orbits around the fixed point. In 
the ideal case, the system would be exactly at the stationary 
state throughout the entire interaction, i.e. 5Pi = 0. We thus 
set the nonlinear adiabaticity condition to be 



'_^{\qi\' + k2\') 












?2 






4 


2 


|(kiP + |^2p) 





5Pi 



< 1 



(25) 



The physical interpretation of r,,/ is as follows. Each of the 
two terms in square brackets represents photon flux excess of 
the low frequency waves over the high frequency waves. The 
first of these terms is for the state under consideration, while 
the second is for the stationary state. Therefore, the complete 
numerator represents the difference in photon flux excess be- 
tween a given set of waves and the stationary state. The de- 
nominator normalizes this quantity by the overall photon flux 
balance between the three waves. 

Using the approximate solution of Eq. |24]for 5 Pi, this con- 
dition becomes 



d{p^/p,) 




d{PflPi) dAT 






dAT d^ 



< V 



(26) 



which means that in order to maintain adiabaticity, the rate of 
change of the normalized stationary state photon flux excess 
in the low frequencies over the high frequencies, P^ /Pt„ has 
to be much slower than the frequency of periodic orbit around 
the fixed point, as expected from classical mechanical theory. 
Eq. |26]is the main result of this work. For the special case of 

l^il 



\q2\ and AF = 0, this inequality leads to 



1 



27^3 



dAT 



d^ 



< 1 



(27) 



Adiabaticity can thus be more closely satisfied when the over- 
all intensity is higher (which increases the overall photon flux 
P3) and when the rate of change of the phase-mismatch is 
lower 

Having established that the system can adiabatically follow 
changes in the phase-mismatch AF, we consider the special 



case where the system is prepared in a nontrivial stationary 
state of Eq. |22l |AF| :s> y/2P^ at the beginning and end of 
the interaction, and AF ends with a sign opposite to the one 
it started with. Clearly, from Fig. [T] when |<7i|^ = |^2|^ the 
adiabatic interaction would result in a complete energy trans- 
fer from the two lower frequencies, COi and CO2, to the high 
frequency, 03. Since it was established that P^ (AF/y^) = 
AF/v^), we will concentrate on adiabatic following 
, where it is readily understood that everything applies 



ofP, 



to P^ upon reversal of the chirp direction. 

In order to demonstrate adiabatic evolution, Eq. [T] were 
solved numerically for three different cases. The results are 
displayed in Fig. |5] In this figure, the dashed curves corre- 
spond to the minus stationary state, calculated using Eq. |22] 
r„i in (c), (f) and (i) was calculated using Eq. |26] In all three 
cases the system started in the minus state. In each case the 
phase-mismatch chirp rate was different, i.e. AF was always 
linearly chirped from — lOy^ to lOy^, but the interaction 
length was varied. In the first case, shown in Fig. H^-c, the 
normalized interaction length was A^ ^/P^ — 1 . Clearly in this 
case the system does not follow the stationary state. Cor- 
respondingly, the adiabatic condition is not satisfied, as r„i 
reaches a value much greater than 1 . In the second case, dis- 
played in Fig. HJl-f, A^ — 10- I" this case the stationary 
state is more closely followed, yet only to a limited extent. 
This is also reflected in the fact that r„i reaches 0.85. Note 
that the area of departure from the stationary state in Fig. |5]l 
and e corresponds to the area where r increases toward 0.85 
in Fig. |5f. Finafly, in the third case, A^^/P^ = 100. Fig. |5^-i 
show that in this case the the stationary state is very closely 
followed, and r„/ ^ 1 throughout the entire interaction. 

In the general case of \qi ^ 7^ 1^2 i-S- P2 7^ 0, Pi will go 
from P3 to 2 |f2 1 ^ ^"3 for increasing AF, when the adiabaticity 
condition is met. This means that energy will be transferred 
from the two low frequencies ©1 and ©2 to the high frequency 
03, until one of the two low frequencies is depleted. A nu- 
merical simulation of such a case is displayed in Fig. |6] where 
P2 = O.3P3. As seen in Fig. |6^, energy is adiabatically trans- 
ferred from the low frequencies to the high frequency until 
none is left at 0)2. From that point on, the three waves in- 
tensities remain essentially unchanged. Fig. ^ shows the 
corresponding value of P1/P3, which indeed goes from 1 to 
(2 IP2I - ^"3) /P3 = -0.4, as expected. 

A special case of the nonlinear adiabatic evolution is the 
case of constant pump approximation'^'^, where the dynam- 
ics becomes linear In this scenario, one of the three waves 
(the pump wave) was taken to be much more intense than the 
other two waves, while another wave was assumed to start 
with no energy. Under the assumption that the effect of the 
interaction on the pump wave is negligible, the remaining two 
waves form a linear dynamical system, to which the linear adi- 
abatic theorem applies. As a result, energy would flow from 
one interacting wave to to other Such a situation was sim- 
ulated here as well, without making the fixed pump approx- 
imation, with the results displayed in Fig. |7] In this case, 
the input pump-to-signal ratio was \q2 (0)|^/ \qi (0)|^ = 100 
and 1^3 (0)1 — 0. Fig. |7^ shows that all of the photon flux 
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Fig. 5. Numerical solutions of Eq. [T]with \qi\^ — |^2p- AF 
is linearly chirped from —IQ^/Pj to lOy^- The system al- 
ways starts in the minus stationary state. The normalized in- 
teraction length is (a-c) A^y/l^ ^ 1 (d-f) A^V^ = 10 (g-i) 
A^v^ = 100. The dashed curves correspond to the minus 
stationary state calculated using Eq. |22] The nonlinear adia- 
batic condition r,,/ in (c), (f) and (i) was calculated using Eq. 
l26l Only the bottom row, in which r„i ^ 1, satisfies the adia- 
batic condition. 
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Fig. 6. Numerical solution of Eq. [T|with the same parame- 
ters as in Fig. |5^, except that P2 — O.3P3. The dashed curves 
correspond to the minus stationary state. 



was transferred from (Oi to 03, with equal contribution from 
(O2 as evident from the inset. This corresponds completely to 
the above description, i.e. the adiabatic interaction took place 
until the coi wave was depleted. Fig. [TJ? shows that P1/P3 
traveled from 1 to 2 |f2| — ^"3 = 0.96, as expected. 

Finally, we note that a trivial stationary state does not corre- 
spond to an elliptic fixed point in the {Qi,Pi) phase space (see 
appendix A), so it would not perform adiabatic following due 
to changing phase-mismatch. This of course can be expected 
on physical grounds, as we do not expect the intensity of the 
only present frequency to be affected by changes in phase- 
mismatch between it and absent frequencies. Interestingly, for 
the case where \qi \^ = \q2\^, each of the two nontrivial station- 
ary states can actually follow the adiabatically-varying phase- 



mismatch into a trivial stationary state with \qi \ 
as evident from Fig. [T]and|4] 
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Fig. 7. Numerical solution of Eq. [T|with the same parameters 
as in Fig. |5^, but assuming a strong pump at 02, a weak signal 
at (Bi and no input energy at 0)3, i.e. the approximate linear 
dynamics regime. The dashed curves correspond to the minus 
stationary state. The inset of (a) shows \q2^ /P-i and has the 
same horizontal axis. 



To summarize this section, adiabatic following can be ob- 
tained when the system is prepared to be near a nontrivial sta- 
tionary state, i.e. such that 8P\ ^ P3, and the rate of change 
of the scaled phase mismatch AF is sufficiently small for the 
given overall photon flux balance P3, as prescribed by Eq. |26] 
If AF changes monotonically, changing signs from beginning 
to end, and |AF| ^ ^/P^ at the beginning and end of the in- 
teraction, the system will evolve adiabatically from Pi = P3 
to Pi = 2 IP2I — P3, or vice versa. The former corresponds to 
upconversion, which ends when one of the two low frequency 
waves is depleted (the one that started with the lower photon 
flux). The latter corresponds to downconversion, which con- 
tinues until the high frequency wave is depleted. In the special 
case where P2 = 0, the system can only evolve from Pi = P3 
to Pi ~ — P3, but not in the reverse direction, since Pi = — P3 is 
a trivial stationary state that does not correspond to an elliptic 
fixed point in the (2i ,Pi ) reduced phase space. 

As a final note, we would like to suggest that the same 
method can be applied to frequency-cascaded and spatially- 
simultaneous TWM processes or higher-order nonlinear adi- 
abatic processes. For example, four wave mixing has also 
been put into canonical Hamiltonian structure, and symme- 
tries, corresponding conservation laws and stationary states 
have been identified^". Optical fiber tapering can be used to 
facilitate adiabatic evolution. A detailed analysis will be car- 
ried out elsewhere. 

B. Bandwidth 

Adiabatic TWM processes have numerically been shown to be 
robust against changes in various parameters, e.g. wavelength 
and temperature t-^iiiii^, which are manifested in changes in 
the phase-mismatch. This robustness stems from the fact that 
AF is swept along a large range of values, so a wide range 
of physical conditions can result in AF within the range that 
satisfies the conditions for adiabatic evolution. 

An estimate of the bandwidth will now be given and 
demonstrated. First, we define the conversion efficiency for 
following the minus state with increasing AF, 
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Fig. 8. Numerically calculated conversion efficiency with 
(a) f 2 = and (b) P2 — Q.3P3, and all other parameters 
the same as in Fig. |5^, for various values of the normal- 
ized phase-mismatch at the center of the interaction medium, 
Ar(0) l\/Pl. The chirp rate and interaction length were kept 
constant. The dashed lines indicate the values of AF (0) / ^JPl 
where the estimation yields r\ — ^. 



Under this definition, r\(P\~P-x) = and 
Tj (Pi = IP2I -^"3) = 1- The full width at half maximum 
of T] is estimated to be (see appendix B for details) 



AFw = AF (A<^ II) - AF (-A^ /2) 



(29) 



The estimated bandwidth is therefore independent of the in- 
tensities of the interacting waves, as it depends only on the 
chirp range of AF. 

The conversion efficiency r\ for f 2 = and Pi = O.3P3 is 
depicted in Fig. [8^ and [8J), respectively, vs. the normal- 
ized phase mismatch at the center of the interaction medium. 
In this simulation, the chirp rate and interaction length were 
kept constant. The vertical dashed lines indicate the locations 
where the estimated efficiency is |, established by introducing 
p- = \P2\ into Eq. [n] For f 2 = and P2 = 0.3^3, the simu- 
lated bandwidth AFbw/V^ is 19.46 and 19.7, respectively.. 
For both cases, the estimated bandwidth is AF^w/V^ = 20, 
which is within 3% of the numerical results. 

For a given chirped phase-mismatch, the bandwidth will 
depend on intensity where intensity determines whether the 
adiabatic evolution conditions are satisfied. On the one hand, 
when the intensity is too low to satisfy the adiabatic condi- 
tion of Eq. |26] the efficiency will always be low. Shifting of 
AF (0) from ^ will more quickly deteriorate efficiency than 
when adiabatic following takes place, thus the bandwidth is 
expected to be lower On the other hand, when the intensity is 
high enough, |AF| ^ y^will never be satisfied, soPj will not 
be close to P3 at the beginning of the interaction. However, in 
this case adiabatic following is still maintained to some extent, 
i.e. the motion of P^ is still slow enough to satisfy Eq. |26] 
so Pi can follow it. Pi will thus orbit the adiabatically mov- 
ing fixed point with a large orbit diameter This will cause the 
efficiency to oscillate rapidly for various AF(0), so a useful 
definition of bandwidth is difficult to find. These phenomena 
are demonstrated numerically in section |4] 

Finally we note that the bandwidth estimation of Eq. [29] 
is valid not only for following the minus state, but when- 
ever the requirements of adiabatic following are satisfied, i.e. 



Pi w f 3 or f 1 w 2 |f2| — ^"3 at the beginning of the interaction, 
AF chirped such that it changes sign from beginning to end, 
|AF| 3> \/Pj at the beginning and end of the interaction and 
Eq. |26]is satisfied throughout the entire process (the details 
can be found in appendix B). It follows that the rest of the dis- 
cussion, regarding intensity too low or too high to satisfy all 
of the aforementioned requirements, is also true for all cases, 
not just those related to P^ and increasing AF. 

4. Numerical Simulations 

In this section, the results of numerical simulations of Eq. [T] 
will be shown, with physical dimensions rather than normal- 
ized units. It will be demonstrated that fully nonlinear, effi- 
cient and wideband adiabatic frequency conversion can read- 
ily be applied in a wide variety of physically available configu- 
rations, using QPM. In all of the simulations presented below, 
the nonlinear medium was taken to be a 40mm long LiNbOj, 
crystal with = 50pm/V^^. The Sellmeier equations of 
Gayer et al.— were used to account for dispersion. 

SFG is addressed first. In this simulation, A2 = 1064.5nm 
and Al is tuned in the range 1450— 1650«m, which yields 
614 < A3 < 647nm. The input intensities of the two low fre- 
quencies are chosen such that they have the same photon flux 
when Al = 1550«m, and the sum frequency wave at A3 was 
always taken to start with no energy. The simulated crys- 
tal had chirped QPM modulation, with a local period starting 
at 11.52/im and ending at 11.79/im. This correspond to AF 
that goes from —3y/Pi to 3y/Pi for a total input intensity of 
200MW/cm2 when Ai = 1550nm. 

Fig. |9^ shows the intensities of the three waves along the 
crystal when Ai = 1550nm and the total input intensity was 
200MW/cm^. As expected, energy is very efficiently trans- 
ferred from the two low frequencies to the high frequency. 
The photon flux conversion efficiency tj, defined by Eq. |28] 
is 0.93. Fig. |9}3 shows the conversion efficiency as a func- 
tion of input wavelength, for several input intensities. For in- 
put intensities of 2 and 200MW /cm^, the maximum efficiency 
was 0.11 and 0.97, with bandwidths of 54.2 and 55.5«m, re- 
spectively. These results correspond to the analysis given in 
subsection 13 Bl the significant increase of efficiency, and the 
slight increase in bandwidth, with intensity, is related to im- 
provement in the satisfaction of the adiabatic condition of Eq. 
l26l For input intensity of 20000MW/cm^, the efficiency per- 
forms oscillations across the Ai tuning range, as predicted, 
due to the fact that the system point is orbiting the fixed point 
from a relatively large distance. 

SHG can be considered as a special case of SFG with 
= llil^, where, additionally, ©1 — 0>2- A simulation 
was conducted for this case well, where the QPM period was 
chirped from 18.83/im to 19.44/xm, once again correspond- 
ing to AF that goes from —3^/Pi to 3^/Pi for input intensity 
of 200MW/cm^ when Ai = 1550nm . All other parameters 
were the same as before. The outcome is displayed in Fig. [TO] 
showing results similar to the case of SFG with (Oi ^ 02. For 
input intensity of 200MW/cm^ , at Ai = 1550nm the conver- 
sion efficiency was 0.96, and the bandwidth was 42nm. 

Difference frequency generation (DFG) is the case where 
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Fig. 9. SFG simulation results: (a) Intensities of the three 
waves along the crystal for input wavelength Xi = 1550nm 
and input intensity lOOMW/cm^ (b) Conversion efficiency vs. 
Xi for different input intensities, which are indicated in units 
of MW/crn^. 



Fig. 11. DFG simulation results: (a) Intensities of the three 
waves along the crystal for input wavelength X2 ~ 1550«m 
and input intensity 200MW fcm^ (b) Conversion efficiency vs. 
X2 for different input intensities, which are indicated in units 
of MW/crn^. 
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Fig. 10. SHG simulation results: (a) Intensities of the two 
waves along the crystal for input wavelength Ai = 1550«m 
and input intensity 200MW /cm^ (b) Conversion efficiency vs. 
Xi for different input intensities, which are indicated in units 
of MW/crn^. 



Fig. 12. OPA simulation results: (a) Intensities of the three 
waves along the crystal for input wavelength X2 = 1550nm 
and input intensity 200MW/cm^ (b) Conversion efficiency vs. 
X2 for different input intensities, which are indicated in units 
of MW/crn^. 



energy is transferred from the high frequency to the two low 
frequencies. In the DFG simulations A3 = 1064.5nm and A2 
was tuned over 1400 — 1800«m, which generates 2605 < Ai < 
4442n;« (consistent with our convention that (Oi < (O2 < (03). 
The QPM period was chirped from 29.86 to 30.86/Xm, and 
here also AF goes from — to for input intensity 

of 200MW/cm^ when A2 = 1550nm (the other low frequency, 
(Oi, always starts with no energy). All other parameters were 
the same as before. In Fig. [TTh it is seen that energy is ef- 
ficiently transferred from the high frequency to the two low 
frequencies, for the case of A2 = 1550nm and input intensity 
of 200MW/cm^. Note that in this case the system follows the 
plus stationary state (Pi starts out negative). The conversion 
efficiency is thus 1 — TJ, which corresponds to the degree of 
depletion of the high frequency. For the case presented in Fig. 
[TTk . the efficiency is 0.95. Fig. [TTb displays the conversion 
efficiency vs. A2 for different input intensities, showing the 
same dependence as in the previous cases. For input intensity 
of 200MW/cOT^, the bandwidth was 212nm. 

The special case of DFG where the input intensity of the 
high frequency is much higher than that of the input low fre- 
quency (0)2 in this case), is commonly denoted OPA. In the 
OPA simulation, the QPM structure was designed such that 
AF goes from —l^/P^ to l^/P^ when A2 = 1550«m and the 
input intensity is 4-QQMW/cm^. Also, the 02 input intensity 
was 20 times lower than the 0)3 input intensity. The resulting 



QPM period was chirped over 29.98 - 30.85^m. All other 
parameters are the same as for the DFG simulation. Fig. [T2t 
shows the intensities of the three waves along the crystal for 
A2 = 1550nm and input intensity of 4-QOMW / cm^ . As before, 
energy is seen to efficiently transfer from the high frequency 
to the two low frequencies, resulting in conversion efficiency 
of 0.97. From start to end, the intensity at 02 was amplified 
by a factor of 13.6. In Fig. [TTb the conversion efficiency is 
plotted vs. A2 for different input intensities, with the same be- 
havior as noted above. A detailed numerical investigation of 
adiabatic OPA has been conducted by Phillips et al.'^. 

The range of parameters used above shows that adiabatic 
TWM can readily be used with nanosecond to picosecond 
pulsed lasers in bulk media or continuous-wave lasers in 
guided structures (e.g. QPM waveguides^). Shorter pulses 
could be stretched, converted and compressed againi"— . 

The combination of broad bandwidth and intensity depen- 
dence of efficiency suggests another application of fully non- 
linear TWM - cleaning the unwanted pedestal of intense ultra- 
short pulses^^'^'*. This could be performed using two QPM 
crystals, as follows. First, the input beam should be linearly 
polarized at 45 degrees to two of the crystals optical axes, 
namely the ordinary and extraordinary axes. In this manner, 
half of the input energy would be at the ordinary polarization, 
and half at the extraordinary polarization. We denote these 
frequency and polarization components o„ and cOe, respec- 
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lively. The first crystal will perform cross-polarized adiabatic 
SHG of the extraordinary polarization, i.e. (Oe + (Oe ^ ^cOo- 
Since conversion efficiency depends on intensity, the high- 
power parts of the pulse will be more efficiently converted 
than the low-power parts. Therefore, after the first crystal, the 
COe wave contains the remaining low-power parts of the pulse. 
These are eliminated by placing a polarizer, aligned along the 
ordinary axis, following the first crystal. After the polarizer, 
we are left with the generated 2cOo wave and the original (un- 
cleaned) co,, wave. These waves now enter the second crystal, 
which performs the degenerate cross-polarized DFG process 
2cOo ~ (Oo ^ (Oe- Once again, the process favors the high- 
power parts of the pulse at (Oo- Passing the output through 
a polarizer aligned along the extraordinary wave will elimi- 
nate the residual low-power at cOu as well as 2cOo, leaving only 
the cleaned pulse at (Oe. 

5. Conclusion 

Adiabatic TWM with fully nonlinear dynamics was put on a 
firm physical basis by rigorous analysis, detailing the condi- 
tions for obtaining adiabatic evolution. Just as the adiabatic 
TWM in the linear dynamics regime was developed from an 
analogy with linear quantum systems'*-^, the method used here 
also follows, in general terms, an analysis of adiabatic evolu- 
tion of nonlinear quantum systems^i*^^*^. Furthermore, the 
nonlinear adiabatic condition was determined, and an estima- 
tion of the bandwidth of adiabatic TWM processes was de- 
rived and shown to be consistent with numerical results. In 
addition, numerical simulations were used to demonstrated 
fully nonlinear adiabatic frequency conversion in several con- 
figurations attainable with current technology. Specifically, 
adiabatic SPG, SHG, DFG and OPA were all shown to be 
efficient over a wide band of input frequencies, using inten- 
sities characteristic of nanosecond pulses in bulk interactions 
or continuous-wave lasers in guided structures. It was also 
explained how adiabatic TWM could be used to facilitate ef- 
ficient pulse cleaning. Finally, it was suggested that adiabatic 
evolution of frequency-cascaded and spatially-simultaneous 
TWM processes or higher order nonlinear processes, such as 
four wave mixing, can also be treated using the same method. 



of motion 



d 



5P, 







(A2) 



Note that the variation in AF causes to be (^dependent, and 
thus dP^ /dt, functions as a source term in Eq. IA2I whereas 
dQ^ /dB, — (see Eq. l2Tl i. We solve for dPi, assuming an 
initial condition of dPi ~ 0, by diagonalizing the coupling 
matrix, which yields 



dp 



cos V I 



where 













dPl 



(A3) 



(A4) 



is the magnitude of each of the two imaginary eigenvalues of 
the coupling matrix ±/v. The nontrivial stationary states fixed 
points are thus elliptic, where v is the frequency of periodic 
orbits around the fixed point. Since it is assumed that the sys- 
tem is near an elliptic fixed point, this frequency is large com- 
pared to all other rates of variation, so the most significant 
contribution to the integral of Eq. IA3l comes from ^ ~ ^' . We 
can therefore approximate 8P\ by taking dPi/dS,' at ^' = ^, 
which can then be taken outside of the integral, yielding 



dp 



1 <iP, 



V d^ 



(A5) 



which was used for Eq. |24] 

Finally, we note that the discussion above referred to a non- 
trivial stationary state. Repeating the same analysis for a triv- 
ial stationary state, for which H — — (AF/S) (Pi + 3P3), yields 
a matrix with two zero eigenvalues in Eq. IA2I Therefore, a 
trivial stationary state does not correspond to an elliptic fixed 
point in the {Qi^Pi) phase space. 



Appendix A: Linearization of the Canonical Hamiltonian 
Dynamics 

This appendix details the linearization procedure that was uti- 
lized to obtain Eq. |24] Linearization of Eq. [TT] i.e. of 
dH/dPi and dH/dQi, can be accomplished in a single step, 
by approximating the Hamiltonian H with a Taylor expansion 



around a fixed point [Qf,P^) up to second order: 



H{Qj,Pj)^H(^Qf,Pf) 



521 + 



dp\ 



(Al) 

where 5Qi =Qi— , 5P\ = P\ —P^, and we have used Eq. 
I20I and also d^H/dQiPi \ = 0. Substituting the ap- 

proximate Hamiltonian in Eq. [T7]leads to the linear equations 



Appendix B: Bandwidth Estimation 

An estimate of the full width at half maximum of the conver- 
sion efficiency tj (see Eq. |28]i will now developed. As noted 
in subsectionEB 77 (fi = ^3) = and T] (Pi = IP2I - ^"3) = 1- 
Furthermore, T] (Pi — IP2I) = 5, i.e. T] — j when Pi is exactly 
half way between P3 and 2 |P2 1 — P3. If Pi starts at P3 and fol- 
lows Pj", it is expected that Pi will end up at |P2|, i.e. half way 
to 2IP2I — P3, if the stationary state fixed point P^ has trav- 
eled the same distance. Assuming a very large chirp range, 
g^uch that Pf always starts near P3 or ends near 2 |P2 1 — P3 (or 
both), there are two cases in which this may happen: (i) P^ 
starts near P3 and ends up at IP2I (ii) Pf starts at IP2I and ends 
near 2 IP2I — P3. In the first case the estimation is more accu- 
rate, since Pi starts near P[^ and will thus follow it as expected 
from the above theory. In the second case. Pi starts near P3 
while Pj" starts at |P2|, i.e. they are not near, so 5Pi ^ P3 
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is not satisfied. Still, as a first order approximation, we can 
expect Pi to traverse a path of similar length to that of . 
Thus, in the first case the condition = jfil is satisfied by 
Ar at the end of the interaction, while at the second case it is 
satisfied at the start. The difference between these two values 
of Ar is, by definition, the chirp range, i.e. the bandwidth is 
estimated to be 

ArBw=Ar(A^/2)-Ar(-A^/2) (Bl) 

The above analysis assumes following of the minus state 
with increasing AF, however it applies to the general case of 
adiabatic following. First, when the minus state is followed 
with decreasing AF, the efficiency is simply 1 — Tj, so the two 
conditions for Tj = | are clearly the same for 1 — rj — j. Fur- 
thermore, following the plus state is the same as following the 
minus state with the opposite chirp direction, so once again 
the same conditions apply. The estimation is therefore valid 
whenever the requirements of adiabatic following are satis- 
fied, i.e. Pi w f 3 or Pi k2\P2\— P3 at the beginning of the in- 
teraction, AF chirped such that it changes sign from beginning 
to end, I AF| ^ the beginning and end of the interaction 

and the condition of Eq. |26] 
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